package com.jl.magic.color;

public class LabUtils {

    // 标准白点 D65
    private static final double Xn = 95.047;
    private static final double Yn = 100.000;
    private static final double Zn = 108.883;

    public static double[] rgbToLab(int r, int g, int b) {
        double[] xyz = rgbToXYZ(r, g, b);
        return xyzToLab(xyz[0], xyz[1], xyz[2]);
    }

    private static double f(double t) {
        return t > 0.008856 ? Math.cbrt(t) : (7.787 * t + 16.0 / 116);
    }

    private static double[] xyzToLab(double x, double y, double z) {
        double fx = f(x / Xn);
        double fy = f(y / Yn);
        double fz = f(z / Zn);

        double L = 116 * fy - 16;
        double a = 500 * (fx - fy);
        double b = 200 * (fy - fz);

        return new double[]{L, a, b};
    }

    private static double[] rgbToXYZ(int r, int g, int b) {
        double[] R = linearize(r / 255.0);
        double[] G = linearize(g / 255.0);
        double[] B = linearize(b / 255.0);

        double X = 0.4124 * R[0] + 0.3576 * G[0] + 0.1805 * B[0];
        double Y = 0.2126 * R[0] + 0.7152 * G[0] + 0.0722 * B[0];
        double Z = 0.0193 * R[0] + 0.1192 * G[0] + 0.9505 * B[0];

        return new double[]{X, Y, Z};
    }

    private static double[] linearize(double v) {
        if (v <= 0.04045) {
            return new double[]{v / 12.92};
        } else {
            return new double[]{Math.pow((v + 0.055) / 1.055, 2.4)};
        }
    }


    // -------------------------
    // CIEDE2000 颜色差异实现
    // -------------------------
    public static double deltaE2000(int[] color1, int[] color2) {
        double[] lab1 = rgbToLab(color1[0], color1[1], color1[2]);
        double[] lab2 = rgbToLab(color2[0], color2[1], color2[2]);

        return deltaE2000(lab1, lab2);
    }

    public static double deltaE2000(double[] lab1, double[] lab2) {
        double L1 = lab1[0], a1 = lab1[1], b1 = lab1[2];
        double L2 = lab2[0], a2 = lab2[1], b2 = lab2[2];

        double deltaL = L1 - L2;
        double deltaA = a1 - a2;
        double deltaB = b1 - b2;

        double c1 = Math.sqrt(a1 * a1 + b1 * b1);
        double c2 = Math.sqrt(a2 * a2 + b2 * b2);
        double deltaC = c1 - c2;

        double avgC = (c1 + c2) / 2;
        double G = 0.5 * (1 - Math.sqrt(Math.pow(avgC, 7) / (Math.pow(avgC, 7) + Math.pow(25, 7))));

        double a1Prime = a1 * (1 + G);
        double a2Prime = a2 * (1 + G);

        double c1Prime = Math.sqrt(a1Prime * a1Prime + b1 * b1);
        double c2Prime = Math.sqrt(a2Prime * a2Prime + b2 * b2);
        double h1Prime = Math.toDegrees(Math.atan2(b1, a1Prime));
        if (h1Prime < 0) h1Prime += 360;
        double h2Prime = Math.toDegrees(Math.atan2(b2, a2Prime));
        if (h2Prime < 0) h2Prime += 360;

        double deltaHPrime;
        if (Math.abs(h1Prime - h2Prime) <= 180) {
            deltaHPrime = h2Prime - h1Prime;
        } else {
            deltaHPrime = (h2Prime - h1Prime + 360) % 360;
            if (deltaHPrime > 180) deltaHPrime -= 360;
        }

        double deltaH = 2 * Math.sqrt(c1Prime * c2Prime) * Math.sin(Math.toRadians(deltaHPrime) / 2);

        double SL = 1.0 + (0.015 * (L1 + L2 - 2 * 50)) / Math.sqrt(20 + Math.pow(L1 + L2 - 2 * 50, 2) / 24);
        double SC = 1.0 + 0.045 * avgC;
        double SH = 1.0 + 0.015 * avgC;

        double RT = Math.sqrt(Math.pow(avgC, 7) / (Math.pow(avgC, 7) + Math.pow(25, 7))) *
                Math.sin(Math.toRadians(30 * Math.exp(-Math.pow((avgC - 27) / 25, 2))));

        return Math.sqrt(
                Math.pow(deltaL / SL, 2) +
                        Math.pow(deltaC / SC, 2) +
                        Math.pow(deltaH / SH, 2) +
                        RT * (deltaC / SC) * (deltaH / SH)
        );
    }

}
